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Abstract 

We consider unsteady undular bores for a pair of coupled equations of Boussinesq- 
type which contain the familiar fully nonlinear dissipationless shallow- water dynam- 
ics and the leading-order fully nonlinear dispersive terms. This system contains one 
horizontal space dimension and time and can be systematically derived from the full 
Euler equations for irrotational flows with a free surface using a standard long-wave 
asymptotic expansion. In this context the system was first derived by Su and Gard- 
ner. It coincides with the one-dimensional flat-bottom reduction of the Green-Naghdi 
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system and, additionally, has recently found a number of fluid dynamics applications 
other than the present context of shallow-water gravity waves. We then use the 
Whitham modulation theory for a one-phase periodic travelling wave to obtain an 
asymptotic analytical description of an undular bore in the Su-Gardner system for a 
full range of "depth" ratios across the bore. The positions of the leading and trailing 
edges of the undular bore and the amplitude of the leading solitary wave of the bore 
are found as functions of this "depth ratio" . The formation of a partial undular bore 
with a rapidly-varying finite-amplitude trailing wave front is predicted for "depth 
ratios" across the bore exceeding 1.43. The analytical results from the modulation 
theory are shown to be in excellent agreement with full numerical solutions for the 
development of an undular bore in the Su-Gardner system. 

I INTRODUCTION 

In shallow water, the transition between two different basic states, each characterized by 
a constant depth and horizontal velocity, is usually referred to as a bore. For sufficiently 
large transitions, the front of the bore is often turbulent, but as noted in the classical 
work of Benjamin and Lighthill IJ, transitions of moderate amplitude are accompanied by 
wave trains without any wave breaking, and are hence called undular bores. Well-known 
examples are the bores on the River Severn in England and the River Dordogne in France. 
Undular bores also arise in other fluid flow contexts; for instance they can occur as internal 
undular bores in the density-stratified waters of the coastal ocean (see, for instance, 0,0), 
and as striking wave-forms with associated cloud formation in the atmospheric boundary 
layer (see, for instance, jl], [Sj). They can also arise in many other physical contexts, and 
in plasma physics for instance, are usually called coUisionless shocks. 

The classical theory of shallow-water undular bores was initiated by Benjamin and 
Lighthill in [Ij. It is based on the analysis of stationary solutions of the Korteweg - 
de Vries (KdV) equation modified by a small viscous term j^i- Subsequent approaches 
to the same problem have been usually based on the Whitham modulation theory (see 
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[7j, [S]), appropriately modified by dissipation; this allows one to study analytically the 
development of an undular bore to a steady state (see 0, [ID], [H])- Most recently, this 
approach was used in ^2] to study the development of an undular bore within the context 
of a bi-directional shallow-water model, based on an integrable version of the well-known 
Boussinesq equations, but modified by a small dissipative term. 

However, a different class of problems arises when one neglects dissipation and consid- 
ers only the dispersion-dominated transition between two basic states. Such a transition 
has the nonlinear oscillatory structure similar to classical undular bores but due to the 
absence of any dissipation, it remains unsteady and expands in time; consequently its qual- 
itative properties and quantitative description are essentially different. This "conservative" 
undular bore evolution can be viewed as the initial stage of the eventual development of 
a steady weakly dissipative undular bore (see, for instance, the numerical simulations of 
unsteady undular bores in a set of Boussinesq equations describing weakly nonlinear long 
water waves in [^, and the analogous simulations for the full Euler equations with a free 
surface in ^1], ^3], ^H])- Importantly, it also represents a universal mechanism for the 
generation of the solitary waves, which can be realized in several different fluid flow con- 
texts. For instance, it occurs in the generation of solitary waves in trans-critical flow over 
topography (see ^Tj, JHl); and in the formation of internal solitary waves in the coastal 
ocean jHj and atmosphere [Sj. If the amplitude of the initial difference between the two 
basic states (e.g. a step in the simplest case) is small enough, the dynamics of an undular 
bore can be described by an integrable equation (typically the KdV equation). In this case 
one can take advantage of the exact methods of integration available for such equations. An 
asymptotic analytical description of an unsteady undular bore was flrst constructed in the 
framework of the KdV equation by Gurevich and Pitaevskii JH] using the Whitham modu- 
lation theory [7j, jH]. The undular bore can then be described by a similarity solution of the 
Whitham modulation equations, obtained by an asymptotic singular reduction from the 
integrable KdV equation. Signiflcantly in this integrable case, the modulation equations 
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can be obtained in the Riemann invariant form. Later, this Gurevich-Pitaevskii theory has 
been generahsed to other integrable equations such as the defocusing nonhnear Schrodinger 
equation [20], [21], [22], the Benjamin-Ono equation [23] and the Kaup-Boussinesq system 
[21]. All these cited works make essential use of the Riemann form of the modulation 
system, a feature associated with the complete integrability of the original equation. 

In recent decades, largely due to extensive theoretical, numerical and experimental 
work on the internal solitary waves (see for instance [23, [2E], |2Z|, |2H|, [20] and references 
therein) it has become clear that although completely integrable systems may successfully 
capture many features of the propagation of weakly nonlinear waves, they may fail to pro- 
vide an accurate description of finite-amplitude dispersive waves. Consequently, significant 
efforts have been directed towards the derivation of relatively simple models enabling the 
quantitative description of the propagation of fully nonlinear waves, and also amenable to 
analytical study. In the context of one-dimensional shallow-water waves such a model was 
derived using a long-wave asymptotic expansion of the full Euler equations for irrotational 
flow by Su and Gardner (SG) [SU]. This system has the same structure as the well-known 
Boussinesq system for long water waves, but importantly for our purposes, retains full 
nonlinearity in the dispersive terms. For convenience, we will present a summary of this 
derivation in an Appendix. Later, it transpired that the SG equations coincide with the 
one-dimensional flat bottom reduction of another model derived by Green and Naghdi (GN) 
[HI] for two-dimensional shallow-water waves. The GN model was obtained using the "di- 
rected fluid sheets" theory, which does not formally require an asymptotic decomposition, 
but instead imposes the condition that the vertical velocity has only a linear dependence 
on the vertical {z) coordinate, and simultaneously assumes that the horizontal velocity is 
independent of z. As a result, the incompressibility condition, the boundary conditions at 
the free surface and at the bottom and an energy equation are satisfied exactly. However, 
as pointed out in [22], this two-dimensional GN model has varying vorticity, even if ini- 
tially the flow is irrotational. Thus, in the context of shallow-water waves propagating on 
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an undisturbed uniform flow, this implication of the two-dimensional GN model contra- 
dicts the exact zero-vorticity requirement of solutions of the full Euler equations for flows 
which are initially irrotational. Thus, the physical validity of the original two-dimensional 
GN model for fully nonlinear shallow-water flows is suspect. This defect in the original 
GN system was recently removed in |33i where the principle of virtual work was com- 
bined with Hamilton's principle expressed in Lagrangian coordinates to derive a new GN 
system (denoted as IGN) valid for irrotational flows. Nevertheless, the one- dimensional 
reduction of the original GN system, being equivalent to the SG system, and also to the 
one-dimensional reduction of the IGN system, does not have this disadvantage. But, al- 
though the SG and GN systems are mathematically equivalent in one spatial dimension, 
their physical meaning is essentially different as the two models are obtained, in fact, by 
quite different approximation methods. 

We also note that the same SG (or GN) system of equations (but for different quantities) 
has appeared in the modelling of nonlinear wave propagation in continua with "memory" 
jHl] , bubbly fluids [321 and in solar magnetohydrodynamics jSH] • With some modifications 
it has also been used for modelling short surface waves governed by the Euler equations 
jin] and for capillary-gravity waves jH]. The "two-layer" generalisation of the same one- 
dimensional fully nonlinear shallow-water system for large amplitude internal wave prop- 
agation has been obtained in ^22i using asymptotic expansions and recently re-derived in 
|HH! using Whitham's Lagrangian approach jSH] • Remarkably, despite the different physical 
origins of the governing equations, the consistent long-wave asymptotic procedure applied 
to this variety of fluid flows yields the same (quite peculiar) form of the nonlinear disper- 
sive terms in the resulting system of equations. This is a clear indication of the universal 
nature of the dispersive term in the SG system of equations. 

Thus, these equations (provided they are used in a consistent asymptotic approxima- 
tion) represent an important mathematical model for understanding general properties of 
large-amplitude dispersive long waves. This is also supported by numerical evidence of 
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the good agreement (superior to that for other generahzed Boussinesq-type models) for 
vertically averaged features of numerical solutions of the full Euler equations in two dimen- 
sions (i.e. one vertical and one horizontal coordinate) describing flow over topography with 
numerical solutions of the GN (i.e. the SG) system (see |i^). 

Note that, although in some of the mathematical literature it is customary to call 
even the one-dimensional version of the fully nonlinear shallow-water equations the Green- 
Naghdi system, here, to avoid confusion with the original, direct Green-Naghdi physical 
model, we will not use this term for the asymptotic system studied in this paper. Instead, 
we will call it the SG system, which would appear to be at least historically more correct. 

While the properties of the permanent steady wave solutions of fully nonlinear systems 
have been intensively studied using analytical means, the properties of unsteady fully non- 
linear dispersive waves and, in particular, the dynamics of dissipationless undular bores, 
remains largely unexplored from an analytical point of view. This absence of "unsteady" 
analytical results, analogous to those which can be obtained in weakly nonlinear wave theo- 
ries, is due to the lack of an integrable structure for most fully nonlinear dispersive systems, 
thus preventing the use of the powerful methods of the inverse scattering theory. In this 
situation, an alternative approach is available through the afore-mentioned Whitham mod- 
ulation theory, which allows one to obtain evolution equations for the local parameters of 
periodic travelling wave solutions. However, unlike the situation for the modulation equa- 
tions associated with completely integrable equations, the modulation equations associated 
with non-integrable wave equations do not possess the Riemann invariant structure. This 
feature is a serious obstacle to their investigation. 

However, recently an analytic approach has been developed in this case ^3], ^1], IIH], 
which allows for the determination of a set of "transition conditions" across an unsteady 
undular bore, and which does not require the presence of the Riemann invariant structure 
for the Whitham modulation equations. This approach takes advantage of some specific 
properties of the Whitham modulation equations connected with their origin as certain 
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averages of exact conservation laws, and so allows for exact reductions in the zero-amplitude 
and zero-wavenumber limits. As a result, the speeds of the boundaries of the undular bore 
region can be calculated in terms of the "depth ratio" across the undular bore, and also 
one then gets the amplitude of the lead solitary wave, a major parameter in observational 
and experimental data. This approach assumes the asymptotic validity of the modulation 
description of the undular bore, an assumption which can be inferred from the rigorous 
results available for completely integrable systems (see jlHI and references therein). It is 
important, therefore, to compare the analytical results of the modulation approach with 
the results of direct numerical simulations in order to establish its validity in non-integrable 
problems. 

In this paper, we use the afore-mentioned Whitham modulation theory to describe an- 
alytically unsteady undular bores in the fully nonlinear shallow water theory. This we will 
do for the whole range of allowed basic states. Then we will make a comparison with the 
full numerical solution of the original system for the same problem. Such a study, while 
elucidating some key properties of fully nonlinear shallow-water undular bores, is also a 
necessary step for understanding the more complicated dynamics of fully nonlinear internal 
undular bores, described by the "multi-layer" generalisation of the system under study in 



II. GOVERNING EQUATIONS AND PERIODIC TRAVEL- 
LING WAVE SOLUTION 

We consider the SG system describing fully nonlinear, unsteady, shallow water waves in 
the form [SO! 



Vt + {w)x = 

Ut + uu^ + rj^ = 



(1) 
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Here rj is the total depth and u is the layer-mean horizontal velocity; all variables are 



non-dimensionalised by their typical values. The first equation is the exact equation for 
conservation of mass and the second equation can be regarded as an approximation to 
the equation for conservation of horizontal momentum. The system (0) has the typical 
structure of the well-known Boussinesq-type systems for shallow water waves, but differs 
from them in retaining full nonlinearity in the leading-order dispersive term. Indeed, 
linearization of the right-hand side of the second equation about the rest state t] = l,u = 
reduces the system ((H) to a familar Boussinesq system (see, for instance, |1Z1)- In Appendix 
A we show how the system (Q) can be consistently obtained from the full irrotational Euler 
equations using an asymptotic expansion in the small dispersion parameter e = ho/ L <^ 1, 
where ho is the (dimensional) equilibrium depth and L is a typical wavelength. We stress, 
however, that there is no limitation on the amplitude. Also we note that the system 
represents a one-layer reduction of the system of Choi & Camassa ^22^ for fully nonlinear 
internal shallow-water waves and also has other physical applications already mentioned in 
the Introduction. The latter can be viewed as an indication that the nonlinear dispersive 
term in the right-hand side of the second equation in has, to a certain degree, a 
universal nature. This suggests that the system ^ is an important mathematical model for 
understanding general properties of fully nonlinear fluid flows beyond the present shallow- 
water application. In particular, although in the context of shallow-water waves the system 
(^, as for any layer- mean model, is unable to reproduce the effects of wave overturning 
and becomes nonphysical for amplitudes greater than some critical value, it is instructive 
to study its solutions for the full range of amplitudes. 

The linear dispersion relation of the SG system (P) relative to a constant basic state 
rio,Uo, has the form, for a right-propagating wave. 



Here c^o is the frequency of the linear waves and k is the wavenumber. 

The weakly nonlinear version of the SG system is obtained by using the standard 




(2) 
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scaling u = 0{6), t] = r]o + 0{6), x = t = 0(5"^/^), where 6 = max[r7 - r]o] ~ 

maxM << 1 is an amplitude parameter. For uni-directional propagation, this leads to the 
KdV equation ior ( = r] — rjo, 

dtC + d.c + ^(d.c + ^C.C = (3) 

(see jSII, m)- 

The first three conservation laws for the system (P) (there exists at least one more - 
see jnij) are obtained by simple algebraic manipulations and represent mass conservation, 
irrotationality and horizontal momentum conservation respectively, 

Vt + ivu)cc = , (4) 



t 



u + -riu^xj + ( y + - y ( + -MMxx - (mx) )) =0, (5) 



{w)t+ + w'^ -^V^ {^xt + uur,^- (u^y)^ =0. (6) 

The periodic travelling wave solution of the SG system (0) is obtained from the ansatz 
1] = rj{6), u = u{9), where 9 = x — ct and c is the phase velocity. Substitution into the 
original system ((T)) and subsequent integration leads to 

where B,C,D are constants of integration. Introducing the roots of the polynomial 

p(^r]) = -{r] - T]i){r] - r]2)ir] - r]3) , Vs > V2 > Vi > ^ , (8) 

we get D"^ = rjirj2rjz- Equation ^ for the total depth r) can be solved in terms of the 
Jacobian elliptic function cn(?/;m): 



r^[e) = r,, + acn'\\\P^^^e-,m\ , (9) 



where ci = V3 ~ V2 , m = — — — (10) 

are the wave amphtude and the modulus respectively. The wavenumber is readily found 
from dZj), © as 

k= ( J^=] =1/^^^^^, (11) 



y/WMU \ J ^/P{r]) / V V1V2V3 2K{m) 



\'?2 



where K{m) is the complete eUiptic integral of the first kind. When m = 1 (772 = Vi) the 
cnoidal wave © becomes a sohtary wave, 



V = VsiO) = ivs - r?i)cosh ^ I V^E^^il 



Next, we define averaging over the periodic family ((7j) by 

F{r]i,ri2,V3,c) = ^ j F{d- r]i, r]2, Vs, c)dd = ^^^^^ j ^^^=y^^ ' ^^^^ 

ri2 

where F{ri) = F{9{ri);rii,ri2,ri3,c). In particular, 

2 

V = , {r]iK{m) + {r]3 - r]i)E{m)) , (14) 

VV3 - Vi 



« = c-^J^=n,(-*^,™|; (15) 



w = vc-\/wmi- (16) 

Here E{m) and Ili{4>] m) are the complete elliptic integrals of the second and third kinds 
respectively. 

Thus, the periodic solution of the SG equations ((1} is characterised by four integrals of 
motion 771, 772, r^a, c or, equivalently, fj, u, k, a. Next, if we allow them to be slowly-varying 
functions of x, t and still require that the travelling wave Q is a solution of the system 
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(to leading order in the small parameter characterizing the ratio of the spatial period 
of the travelling wave Q to the typical scale for its variations in space) we arrive at the 
system of Whitham modulation equations. 

III. WHITHAM MODULATION EQUATIONS 

There are several methods to obtain the modulation equations. We will follow here the 
original Whitham method of averaging the conservation laws [Zj, which is equivalent to 
formal multiple-scale asymptotic expansion but is more convenient for our purposes. 

To this end, we apply the averaging |T3|l to the conservation laws (jH) - (jH} considered 
for the periodic family 0. Due to the characteristic scale separation for the periodic so- 
lution and the modulations, the operations of differentiation and averaging asymptotically 
commute, and so we arrive at a set of quasilinear equations which can be represented in 
the form 

d— _ _ d — _ _ 

—Pj{T],u,k,a) + —Qj{r],u,k,a) = 0, j = 1,2,3. (17) 

The system ()17p should be closed by the wavenumber conservation law, which serves as a 
consistency condition in the formal asymptotic procedure equivalent to the Whitham av- 
eraging method (indeed, this equation can be obtained by averaging an extra conservation 
law and combining it with the averaged equations (fTTjl [7], jS]), 

d d _ _ 

^/j; + ^^(^?, M, fc, a) = . (18) 

The frequency u in (fTHj) is defined as 

u = kc{f],u,k,a) , (19) 

where the phase speed c is expressed in terms of the basic modulation variables with the 
aid of Eqs. (|TO|l . (fTT| . (fTH) . (fT5|l . We recall that typical x, t-scale for variations of the 
dependent variables in ()17j) . ()18|1 is much larger than that for the periodic travelling wave 
(ini) itself in the original equations ((T)). 
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Using (IH) - (0), (fTT|) . (Uni), we can readily obtain explicit expressions for the densities 
Pj, k and fluxes Qj, u of the Whitham modulation system, in terms of the original param- 
eters 7]j, c which are more convenient to use as intermediate modulation variables. These 
expressions involve complete elliptic integrals similar to the KdV case (see |S] for instance) 
but are, as expected, more cumbersome. The dependence of Pj, Qj on f/, u, k, a in (fTTj) is 
then specified parametrically through Eqs. (fTnjl . (|TT|l . (|T^ . (fTH|) . 

Since we will use only some asymptotic properties of the Whitham system in the small 
amplitude and small wavenumber limits, we do not need to present the full expressions 
here. It is sufficient here to note that the modulation system of equations can be re-written 
in a generic quasi-linear form, 

y, + B(y)y,. = 0, (20) 

where y = {f],u,k,a)'^ and the entries of the coefficient matrix are defined by -Bjj(y) = 
dy-Qj/dy.Pj, i,j = 1, ... ,4, and we define P4 = k, Q4 = u . In some special cases, the 
Whitham systems can be represented in diagonal Riemann form, despite the fact that the 
number of the dependent variables exceeds two. This remarkable fact was first established 
by Whitham [2j for the third-order KdV modulation system, and then generalised to many 
other completely integrable systems (see jSO] for a simple method of obtaining the Whitham 
system in Riemann form for a large class of integrable dispersive equations belonging to 
Ablowitz-Kaup-Newell-Segur (AKNS) hierarchy). The presence of the Riemann invariants 
dramatically simplifies further analysis of the modulation equations, and makes readily 
available many important particular solutions (see [21], [22]; [2^ for instance). However, 
for non-integrable equations such a structure is typically not available, which makes the 
corresponding analysis of the modulation system far more complicated. 

There is no indication that the system (0) is integrable so the Whitham system (fT7|) . 
|T8|) is not likely to possess Riemann invariants. One can, in principle, calculate the char- 
acteristic velocities Xj ,j = 1, . . .4 for this system specified by the roots of the determinant 
det(B — AI) = (where I is the unit matrix) but in the absence of the underlying algebraic 
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structure, these expressions are unlikely to be amenable to analytic treatment. Instead, 
we can take advantage of some general properties of Wliitliam systems for obtaining the 
main quantitative features of the solutions, available even in the absence of integrability 
properties. 

We now outline the properties of the Whitham system (fT7j) . (fTHj) that distinguish it 
from the general class of hyperbolic quasilinear systems of fourth order. Most importantly, 
the Whitham system admits exact reductions for the linear a = and solitary wave k = 
regimes. It is clear that, since in both limits the oscillations do no contribute to the 
averaging (|T3|) . one has rj = f], u = u, rfu = fju, and three averaged conservation laws (|T7jl 
must reduce to the dispersionless shallow water equations for fj and u, 

Vt + ivu)^ = , ut + uu^ + r/^. = . (21) 

At the same time, in each of the indicated limits, two out of the four characteristic families 
of the Whitham system must merge into a double characteristic to provide a consistent 
reduction to a hyperbolic system of a lower (third) order. Thus the limits a — and k 
are singular ones for the modulation system. The detailed description of this reduction 
for a general class of nonlinear weakly dispersive systems can be found in Of course, 
for a specific system all the described properties can be established by a direct asymptotic 
analysis of the system (fTTjl . (fT^ . 

In the linear limit (m 0), for c > we get from Eqs. (fTUj) . pi|) . (fT^ . (fT3jl the mod- 
ulation dispersion relation for linear waves propagating on the slowly varying background 
r/(x,t), u{x,t), 

f ^1/2 \ 

a = 0: = uJoiv, u, k) = k + ^ ^^^^^^ j . (22) 

As expected, this is just the linear dispersion relation (j2)) for the right-propagating linear 
wave, where Uq ^ u, tjq f/. 

In the solitary wave limit (m 1), we have from (fT^ . (fTH|) the speed- amplitude relation 
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for a solitary wave propagating about the mean background 77, u\ 



k = 0: 



c = 



Cs{ri, u, a) = u + ^Jr] + a . 



(23) 



We note that for the constant background flow formula (j^Hj) appears in Rayleigh [SI] . It can 
be shown quite generally (see [l^j) that the linear group velocity dojo/dk, and the solitary 
wave speed coincide with the multiple characteristic velocities of the full modulation 
system fITTj) . (fTHj) in the limits a ^ and A; respectively, again, as expected. 

One should emphasise here, that in contrast to the traditional analysis of linear and 
solitary wave propagation, the background flow parameters f/, u in this modulation theory 
are not constant in general, but vary slowly in x, t over a typical modulation scale (Ax ~ 
At ^ 1) according to Eqs. (|2T|l . This x, t-dependence can be easily (parametrically) taken 
into account in the local initial- value problem analysis of the linear wave packet or solitary 
wave train propagation, which makes it essentially equivalent to the constant background 
flow case. However, in the undular bore problem, where the solitary waves and the linear 
wavepacket are restrained by being the parts of a global nonlinear wave structure, this 
dependence plays a crucial role. 

We now present a weakly nonlinear asymptotic reduction of the modulation equations 
(fTTjl . (fT^ . which will be important for our further analysis. One should emphasise that 
the small-amplitude expansion of the modulation system for equations (0) will not lead 
directly to the much-studied KdV modulation dynamics as the KdV asymptotics Q, along 
with small amplitudes, implies also a further long- wave scaling. We expand (fT7j) . (fT8|l for 
m ^ 1 and express the result, using the representation of Eqs. ()13|l - (fT?)|l in terms of 
physical variables, rather than polynomial roots rjj. After a simple but somewhat lengthy 
calculation we obtain for a <^ 1 (see for a similar representation for the KdV modulation 
system and |^ for nonlinear plasma wave modulation equations), 




(24) 
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dt dx\2 



— + — {uoiv, u, k) + u;2{v, k)A^) = 0{A^d,r], A'd.k, A'd.A') . (27) 

Here ujQ{f],u,k) is obtained from ujQ{k) given by Eq. (0) with tjq and Uq replaced with fj 
and u as in (j22|). The values A"^, /i, U2 are expressed in terms of t], k, a by the formulas 

=W-^^= 8^1/2(1 + ^2)1/2 (28) 
^=2r,V2(i + , 2)3/2' ^^2 = ^ 8«:2(l + ^2) ' = ^ ^ • (29) 



The wave energy equation ()26p is obtained by subtracting the modulation equations 

and ()25|l . multiplied by u and respectively, from the averaged momentum conservation 

equation {f]u)t + {. . .)x = 0. It is essential that the terms 0{A'^dxk) do not appear in the 

right-hand part of equation (j^Bj) . i.e. the coefficient for d^k in the considered approximation 

is exactly A'^d'^uJo{f],u, k)/dk'^. 

One can see that the system (j^^ - (j2Z|) indeed admits an exact reduction for A = 

(which is to say a = 0) to the ideal shallow-water equations ()21|) for fj, u complemented by 

the wave conservation equation in the zero- amplitude limit: 

A = 0: dk^ duo{n.u,k) ^ 
at ox 

An analogous asymptotic analysis can be performed for the solitary wave limit, when 
(1 — m) ^ 1 (i.e. k <^ 1) resulting, for the limiting case A; = 0, in the same shallow- 
water reduction complemented, instead of by the equation for the solitary wave 
amplitude, which can be represented in the general form [52] 

da _ da p \du r /- / \ 

k = 0: — + Cs{u,f],a)— + fi{f],a)— + f2{r],a)— = 0. (31) 
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where /i, 2(^,0') are certain functions which can be, in principle, obtained by passage to 
the singular limit /c ^ in the modulation system ()17p. pSj). The actual calculation of 
this limit is significantly more cumbersome than in the linear case as one has to make 
the approximations up to exponentially small terms ~ exp(— l//c). We shall not derive 
the amplitude equation (jHT|) here explicitly; instead, we will later take advantage of an 
alternative (conjugate) system of modulation variables enabling one to perform a singular 
zero-wavenumber limiting transition directly in the wavenumber conservation law (|T^ . 

IV. UNSTEADY UNDULAR BORE TRANSITION 
A. General construction 

We shall consider initial conditions for the system (Q) in the form of a step for the 
variables rj and u: 

t = : T] = T]^ , M = M~ for X < 0; t] = 77"^, m = m"*" for x > , (32) 

where 77^ and are constants. Without loss of generality, we will later set rj^ = 1 and 
u'^ = for convenience. 

We shall model the large-time asymptotic structure of the undular bore with the aid 
of an expansion fan solution of the Whitham system (fT7j) . p8|l for the local parameters 
of the travelling wave ((Zj). This approach to the modulation dynamics was first proposed 
by Gurevich and Pitaevskii for the KdV equation, and has proved to be very effective 
for description of weakly nonlinear unsteady undular bores (see for instance ^\ , ^H] , and 
[3j). The Gurevich-Pitaevskii formulation has been extended to other integrable models 
such as the defocusing nonlinear Schrodinger equation j20|, 10]) > Benjamin- 
Ono equation j^HI and the Kaup-Boussinesq system [21], jSH]- The original Gurevich- 
Pitaevskii formulation makes essential use of the availability of the Riemann invariants for 
the modulation system. Its generalisation to the case when the Riemann invariants are not 
available has been recently developed in [l^, 03] and it is this approach that will be 
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used in the present work. 

We shall first assume that the modulation equations are hyperbolic and genuinely non- 
linear for the solutions of interest, which implies that: 

a) the characteristic velocities of the modulation system (fTTj) . (fTHj) , Ai > A2 > A3 > A4 are 
real and distinct for the solution of interest everywhere except those special lines in the 
X, t-plane where two of them collapse into a double eigenvalue of the Whitham system; 
a) the characteristic fields dx/dt = Xj{f],u,k,a) are not linearly degenerate, i.e. 
Tj ■ gradAj 7^ 0, where rj{fi,u,k,a) is the right eigenvector of the coefficient matrix B in 
(f20|l corresponding to the eigenvalue Xj [5lj . 

The first condition ensures the modulational stability of the solution under study. Two 
properties of the SG system (P) relevant to this issue can be readily inferred from its gen- 
eral structure: 

i) the linearised eigenmode problem for the system (P) yields pure real eigenvalues (the 
frequencies in the linear dispersion relation Q are real for all values of k); ii) the disper- 
sionless limit of full equations is a strictly hyperbolic shallow-water system. 
Also we note the recent result in where the stability of solitary waves with respect 
to linear perturbations has been established using the Evans function technique for a one- 
dimensional Green-Naghdi system mathematically equivalent to (0). All this, however, does 
not guarantee the modulational stability of nonlinear travelling waves with < m < 1. 
This can be established by the analysis of the characteristic velocities of the nonlinear 
modulation system ()17p for the class of solutions of interest. The general condition for 
hyperbolicity of the modulation system for the one- dimensional Green-Naghdi equations 
in Lagrangian variables has been derived in [HI], but the actual verification of this con- 
dition represents an involved technical task, possibly not resolvable by analytical means. 
Later we will establish and explicitly verify a more particular criterion of the modulational 
stability of the undular bores in the system ((T)) using the asymptotic system - (j2Z|). 
As to the genuine nonlinearity issue, we will show below that the modulation system for 
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Eqs. admits linear degeneration of one of the characteristic fields for a certain domain 
of the initial conditions ()32|). and will derive an important restriction connected with the 
formation of partial undular bores in Section 5. 

According to the general Gurevich-Pitaevskii scheme we replace the original initial- 
value problem ()32|) for the system with a matching problem for the mean flow at the 
free boundaries (i.e. the edges of the undular bore). Speciflcally, we require the continuity 
of fj,u along the lines x = x~{t), x = x^{t) deflned by the conditions of zero amplitude 
(traihng edge) and zero wavenumber (leading edge) respectively, i.e. 

x = x^(t): a = 0, f) = ri~ , u = u~ , 

^ ' (33) 

x = x^{t): k = 0, fj = Tj^ , u = u^. 

The edges of the undular bore at a; = x^{t) represent free boundaries, and their deter- 
mination is a part of the solution to be obtained. Since for the step-like initial data ()32|1 
the problem reformulation in terms of the averages does not contain any parameters of the 
length dimension the solution of the quasi-linear Whitham equations ()2m) must depend on 
the similarity variable s = x/t alone. The boundaries x^{t) then represent straight lines 
X = s^t, where are the speeds of the undular bore edges which should be found in terms 
of the given initial step parameters. 

The similarity solution of the Whitham equations can be represented in the general 
form 

Fi{f],u,k,a) = li, Xj{f],u,k,a) = s , i = 1,2,3, (34) 

where are constants and Xj is one of the characteristic velocities of the modulation system 
(fTTjl . ((IHI). They must be chosen so that the solution satisfles the matching conditions 
As a result, the solution (jH^ yields a slow x, t-dependence of the parameters fj,u,k,a of 
the periodic solution Q so that after substitution of into Q we arrive at a slowly 
modulated nonlinear wave degenerating into a linear wavepacket near the trailing edge, 
and gradually transforming into a chain of successive solitary waves close to the leading 
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edge (see figure 1). This wave structure corresponds to the observable features of unsteady 
undular bores (see, for instance, ISE] and jS] for laboratory and natural observations and 
[T!^ . [Ti] . [TK] . pT)] for numerical simulations), and has been rigorously recovered as the 
small-dispersion asymptotics of the exact solution of certain integrable equations (see the 
review |^ and references therein). 



IT 




Figure 1: Oscillatory structure of the undular bore evolving from an initial step (dashed 
line) . 

Explicit modulation solutions in the form of an expansion fan (j3l|l have been con- 
structed for a number of integrable (by the inverse scattering transform) bi-directional 
equations where the functions Fi represent Riemann invariants of the modulation system 
(see for instance [23 and references therein). We emphasise, however, that the general 
similarity solution ()34|1 can, in principle, be constructed for any quasi-linear hyperbolic 
system and, as such, is available for the Whitham system (pUj) regardless of the existence 
of the Riemann invariants. Therefore, the Gurevich-Pitaevskii construction in fact does 
not rely on the integrability of the system in the sense of the availability of an inverse 
scattering transform. 

The fact that the modulation system admits exact reductions as a ^ and — + to 
the dispersionless limit equations ()21|1 . complemented by the wave number conservation 
law in the corresponding limit (see Section 3) allows one to study the limiting structure 
of the integrals in the expansion fan solution ()34|) without necessarily constructing the 

19 



solution itself. Remarkably, these limiting integrals can be expressed in terms of the linear 
dispersion relation of the original system, and one of its "dispersionless" nonlinear char- 
acteristic velocities. The corresponding analysis has been recently described in 
[13] for a rather broad class of equations that one may characterize as 2 x 2 strictly hyper- 
bolic modified by weak dispersion. As a result, a set of the transition conditions has been 
obtained, which allow one to fit the undular bore into the classical solution of the shallow 
water equations much as the classical shock is fitted into the ideal gas dynamics solution 
without constructing the detailed solution in the transition region. Next we derive explicit 
expressions for these conditions for the undular bore in the SG system (P). Our analytical 
results are compared with numerical simulations of the full system (^. An outline of the 
numerical method is given in the Appendix B. 

B. Simple undular bore transition relation 

Since the number of the constants in the solution is three, whereas the initial con- 
ditions (fH^ are characterised by four constants, there should be an additional relationship 

$(r/+,r/-,M+,u-) = 0, (35) 

for the admissible values of the variables at the edges of the undular bore. This relationship 
has been shown in |33] , [IHl to follow from the continuous matching of the characteristics of 
the Whitham system ()17j) . (fT^ and the ideal shallow water system (12111 for the "external" 
flow at the undular bore boundaries. Such a reformulation of the problem in terms of char- 
acteristics is equivalent to the original matching conditions p3|l and yields the relationship 
flH5|l in the form of the zero-jump condition across the undular bore for one of the classical 
Riemann invariants. For the right-propagating shallow-water undular bore it assumes the 
form 




20 



Generally, the step initial condition ()32|) resolves into a combination of two waves: 
an undular bore(s) and/or a rarefaction wave(s). The family of the initial steps which 
resolve into a single right-propagating undular bore described by the expansion fan solution 
of the Whitham equations is distinguished by the condition (jSBj) . which coincides with 
the relationship between u and r] for the simple wave of compression in the limit of the 
dispersionless equations (see [Hj for instance). We shall call the undular bores satisfying 
condition pHjl . simple undular bores. Of course, for a left-propagating simple undular bore, 
the transition relation is obtained from the zero jump condition for another analogous 
classical shallow- water Riemann invariant. 

Without loss of generality one can put in m"*" = , 77+ = 1 so that the simple 
undular bore transition curve becomes 

u- = 2(^/^-1) (37) 

and yields all the admissible states upstream of the undular bore provided the flow down- 
stream is fixed as above. We emphasise that this simple undular bore transition condition 

1.8 
1.6 
1.4 
1.2 
1 

0.8 
0.6 
0.4 
0.2 



Figure 2: Riemann invariant (dashed line) and the shock-wave (dotted line) transition 
curves 

does not coincide with the classical jump relation obtained from a combined consideration 
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of the balance of mass and momentum across the bore provided its width is constant (which 
is the case both for turbulent bores and for established frictional undular bores - see for 
instance [H!;[I2]): 



The curves ()37p and ()38p are shown in figure 2 and demonstrate high contact for small 
jumps, which is expected in view of the well-known fact that the Riemann invariant has 
only a third order jump across a weak shock. Thus, the distinction between the jump 
condition ()38|) and the simple undular bore transition relation becomes noticeable 
only for large-amplitude undular bores. This distinction may not seem very important 
from a practical point of view as the shallow-water undular bores are known to exist only 
for 7]^ /f]^ ^ 1.3 (see [HI, ^M) after which they become turbulent with wave breaking. 
However, it is the simple undular bore curve (rather than the jump condition (jHHj) ) 
that is consistent with the modulation system, and allows analytic determination of the 
undular bore location. Also, equations appear in other physical systems as well as in 
the shallow water context (see [33], [33], [SHj), where different physical restrictions may 
apply and a wider range of amplitudes may be involved. 

C. Undular bore location and lead solitary wave amplitude 

Let this similarity solution of the Whitham equations be confined to the interval s~t < 
X < s~^t so that the modulation provides a transition of the underlying travelling wave ([HJ) 
from a linear wavepacket (a = 0) at the trailing edge to a solitary wave {k = 0) at the 
leading edge. Then the speeds of the undular bore edges can be determined by the 
general expressions derived in [13], ^3], [12] • Below we briefiy outline how this method 
applies to the SG system ([T]). 

The analysis in [13], [H], [12] is based on the general fact that, due to the degeneration 
of the three averaged "hydrodynamic" conservation laws ()17|) for a = and k = into 
the second-order dispersionless system (the shallow-water system ()2H1 for f/ and u in our 




(38) 
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case) the boundaries of the undular bore necessarily correspond to multiple characteristics 
of the Whitham modulation equations. The corresponding multiple characteristic speeds 
can then be obtained without a full integration of the Whitham system by deriving the 
solution in the form (jS^). To this end, one complements the matching conditions (jHHI) 
by the definition of the edges using natural kinematic conditions at the free boundaries 
X = s^t. Namely, at the trailing edge we require that the edge speed to coincide with the 
linear group velocity, and at the leading edge it coincides with the leading solitary wave 
speed, that is 

s~ = ^{ri',u',k') , s+ = c^(r7+,u+,a+) , (39) 

where ujo^fj, u, k), Csifj, u, a) are given by (j221), (|2SI)- These kinematic conditions contain 
unknown parameters k~ and a+, which need to be found in terms of the given initial jumps 
for rj and u. 

The described degeneration of the Whitham modulation equations in the limits a = 
and k = implies the existence of two families of integrals: {a = 0,^i{u,fi) = Ci, 
^2{k,f]) = C2} and {k = 0,^l{u,f]) = Q, $2(0,77) = where Ci^2, C'1,2 are constants. 
In the context of the similarity solutions these integrals represent the zero-amplitude and 
zero-wavenumber "sections" of the general integrals Fj in the full solution ()34|). However, 
the "limiting" integrals $1,2, '^'i 2 can be found directly from the reductions of the Whitham 
equations as a = and A; = 0. The constants C1.2, CI 2 are then found from the matching 
conditions and the transition curve ()3(j|l . As a result, we find that the parameters k~ 
and a"*" are obtained from the functions k{fj) and a(f/) at f/ = 77" and f] = rj^ respectively. 

Indeed, for a = the Whitham system (|T7j) . ()18|1 reduces to a third-order system 
consisting of the shallow-water equations for the mean flow ()21|) complemented by the wave 
conservation equation (jHUj) for linear waves. Then, substituting u = u{ff) and k = k{ff) 
into (PT|) . (jHUj) we obtain the integrals 

df] u{ff) ^ ^ — dujQ / dk 
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It is shown in that the ordinary differential equation in pUj) is just the equation for 
the characteristic integral of the full modulation system along the multiple characteristic 
on which a = and thus, its integral specifies, from the first equation in Eq. (jHO)), the 
speed of the trailing edge of the undular bore. The signs and constants of integration in 
PUj) are then found from the transition relation ()36|) and the matching conditions ((221) ■ 

The reduction of the Whitham system for k = consists of the shallow-water system 
(PT|) and the amplitude equation (j!TT|) . Analogously to the zero-amplitude case, one can 
obtain the integrals u{f]), a{ff) and then, from the second equation in Eq. (jH^ . the leading 
edge speed of the undular bore. There is, however, the afore-mentioned technical compli- 
cation connected with obtaining the amplitude modulation equation (j!?T|l from Eqs. (fTTj) . 
(fTH|) in the limit 0. This difficulty is bypassed in the papers of 112], by the 

introduction of the conjugate wavenumber 

/ m \ -1 

^ = ^1/^1 (41) 




instead of the amplitude a in the modulation equations and considering the (singular) 
limiting transition as /c ^ in the wavenumber conservation law (fTH|) . Then, the equation 
for the characteristic integral k{f]) turns out to coincide with the equation for its zero- 
amplitude counterpart k{f]) (jlUj) . but the linear dispersion relation u!o{f],u,k) should be 
replaced with the conjugate expression ojg = —iuQ{f],u,ik). 

We now apply this construction to the description of the undular bore transition in 
fully nonlinear shallow- water system (P). Again, for simplicity, we put = 1, u+ = in 
the downstream flow. Then the expressions for the speeds of the edges of the undular bore 
take the form 

--tcr.n, (42, 
Here the functions floiv, k) and VLs{fj, k) can be expressed in terms of the linear dispersion 
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relation for the modulations 

noiv, k) = Uo{f],u{f]),k) , ns{f],k) = -ino{f],ik) , (43) 

where u{f]) = 2{^/f] — 1), that is 

Qoiv, k) = 2k{n'/' - 1) + ^^^fj!,]^y/, . (44) 

The parameters k~ and k~^ in (jl^ are calculated as the boundary values k~ = k{ri~), 
k~^ = k{l) of the functions k{f]) and k{fi), which in turn are defined by the ordinary 
differential equations 

dk dilo I df\ 



df] V{r]) - dQo/dk 
dk dVLs/dfj 



fc(l) = 0, (45) 



k{ri-) = Q, (46) 



d7] V{f]) - dVts/dk 
where 

Vij]) = u{f]) + 77^/2 = 3^1/2 _ 2. (47) 

Substituting (jiijl into (jlKjl and introducing a = (1 + k'^ff /?))^^^'^ as a new variable instead 
of k we obtain the ordinary differential equation, 

dfi 2(1 + a + a^) , , , , , 

— = , — ^-rda, a{l) = l. 48 
f] a(l + a)(a-4) ' ^ ^ ^ ^ 

Since (pSj) has a separated form, we can integrate it to obtain 

^ / A \ 21/10 /, , \ 2/5 



Next, using ()42|1 . we get an implicit expression for the trailing edge s~ in terms of 
the total depth ratio across the undular bore A = 77^/77+ = 77": 
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Similarly, using the conjugate linear dispersion relation (p3j) and the equation H4bii we 
obtain from (jUJ, the equation of the leading edge s"^ = s+(A) in an implicit form, 

A / 3 \''/'V 2 



2/5 



We emphasise that Eqs. (jHUjl . (jHT|) are exact solutions of the modulation equations for a 
fully nonlinear wave regime. 




Figure 3: The leading s+ (upper curve) and the trailing (lower curve) edge speeds vs 
the depth ratio A across the simple undular bore. Dashed line: modulation solution (jnOI), 
()5H1 : Diamonds: values of s~^, s~ extracted from the full numerical solution; Dotted line: 
the KdV modulation solution (jSHl). 

The graphs of s"*", s~ versus the depth ratio A across the bore are shown in figure 3 
(dashed line) and demonstrate excellent agreement with the results of direct numerical 
simulations (solid line) of the system (P) for A < Acr ~ 1.43 where Acr corresponds to the 
minimum of the function s~(A). We will discuss the reason for the discrepancy between 
the modulation predictions for the leading edge and numerical results for large depth ratios 
A > Acr in Section 6. 

Next, we consider the correspondence between these exact expressions (jHTjl . (jHUjl for a 
fully nonlinear undular bore, and the weakly nonlinear KdV asymptotics ()53|) . Thus, we 
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introduce in (jSHl) , (EH) a small parameter 5 = A — 1-Clto obtain the expansions 

8-^ = 1 + 6-^6' + 0{6') , ^" = 1 - + f + 0{6') . (52) 

On the other hand, for the KdV equation we have from the original Gurevich-Pitaevskii 
(1973) solution, taking into account the coefficients in the KdV equation Q, 

4dy = l + (A-l)> W = 1-^(^-1)- (53) 

Thus, the expansions agree to first order with the KdV formulas (jSSl)- 

It can be seen from the numerical solution shown in figure 3 and from the asymptotic 
expansion (jK^ of the analytic solution (jKT|) that, for large initial jumps the leading edge 
of the undular bore in the fully nonlinear model (0) demonstrates a tendency towards a 
slower speed compared with its weakly nonlinear counterpart. This feature agrees with 
the results of numerical simulations of undular bores in other fully nonlinear models of 
Boussinesq type and of the full Euler equations for potential flows [12], [TEl- 

At the same time, one can see that the trailing edge of the undular bore in the fully 
nonlinear shallow-water dynamics described by the SG system is noticeably shifted 
forward compared with that of the KdV solution. As a result, the fully nonlinear undular 
bore transition is significantly narrower than that predicted by the weakly nonlinear theory. 
This feature was noted in (16j, where shallow- water undular bores were studied using a 
fully nonlinear numerical model without further long-wave expansions. This result that the 
SG system (^, while formally derived as a long-wave model, can adequately reproduce the 
features associated with the propagation of short waves (the dynamics of the trailing edge 
of the undular bore), is due to the fact that system (0) has a better approximation of the full 
linear dispersion relation than the KdV equation. Indeed, the linear dispersion relation Q 
can be regarded as a Pade approximation of the exact expression {ujq — ku^Y = kta.nh.kriQ 
following from the Euler equations, whereas the KdV dispersion relation is a less accurate 
Taylor series approximation. Applications of modifications of system (PJ to the description 
of the nonlinear dynamics of short waves have been considered in jlU], 
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Figure 4: The lead solitary wave amplitude vs depth ratio A. Solid line: numerical solution; 
dashed line: modulation solution (j54j) 

From the expression ()23|1 relating the speed and the amplitude of the SG solitary wave, 
and the definition of the leading edge s+ = Cs{ri~^, a~^) we derive the equation for the 
amplitude of the leading solitary wave for the free surface elevation, 



The expansion of (jK^ for small jumps 6 = A — 1 yields a+ = 26 + 6'^/6 + 0{6^) which agrees 
to first order with the well-known KdV expression a^^y = 2(A — 1) obtained by Gurevich 
and Pitaevskii JHl- The graph of a"*" versus the depth ratio A is shown in figure 4 (dashed 
line) and agrees with the results of direct numerical simulations (solid line). Again, the 
discrepancy for A > Acr ~ 1-43 will be discussed in the next section. 

One should bear in mind that formulae (jSUj) . ()51|). and ()54|) were obtained for simple 
undular bores i.e. they are valid only if the initial levels in t] and u satisfy the Riemann 
invariant condition (jHTj) . For an arbitrary initial discontinuity the resolution will occur via a 
combination of two simple waves (undular bore(s) and/or rarefaction wave(s)) propagating 
in opposite directions in a certain reference frame. The relevant combination is found by 
determining the pair of intersecting transition (simple-wave) curves in the rj — u plane with 




(54) 
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their centres at the initial basic states, in a similar manner to the construction of decay 
diagrams in gas dynamics (where one of the waves could be a classical shock satisfying the 
Rankine-Hugoniot curve); see [HZj for instance. An important example of such a two-wave 
resolution occurs in the dam-break problem and will be considered in Section VII. 

The construction of the transition conditions for a simple undular bore, described above, 
is subject to the inequalities 

— v^rp < < + \/ri^ , + -\/?/+ < s'^ , > . (55) 

These inequalities are analogous to the entropy conditions of classical gas dynamics (see for 
instance ^^) and ensure that only three of the four "dispersionless" characteristics families 
of the shallow water system 1)211) : x/t = ± Vv^ transfer initial data ()32)1 from the x-axis 
into the undular bore domain in the (x, t)-plane providing consistency with the number of 
parameters in the similarity solution ()34|) (see 145^ for details). A direct verification shows 
that conditions are satisfied for all values of A. 
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Figure 5: The numerical solution of the system for the initial depth ratio A = 1.2, 
t = 250. Initial step data: r/" = A = 1.2, = 2{^/rr - 1), r/+ = 1, m+ = 

Next, we briefly describe the results of the direct numerical simulations used to obtain 
the graphs in figures 3,4. The numerical method itself is outlined in the Appendix B. 
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A typical numerical solution for an undular bore is shown in figure 5, where we plot 
the surface elevation 77. The initial profiles for 77 and u are chosen so that the simple 
undular bore condition (jHTj) is satisfied for the asymptotic states at infinity. As a result, 
a backward propagating rarefaction wave does not form, and the resolution of the step 
occurs through the undular bore alone. In contrast with the model pattern shown in 
figure 1 the numerically obtained solution contains small-amplitude oscillations behind the 
undular bore. These are essentially linear oscillations arising due to the deviation of the 
initial perturbation used in the numerical simulations from the Heaviside step function. 
Indeed, these oscillations decrease in amplitude proportionally to t^^/^, as expected in a 
linear theory. A detailed analysis of these linear oscillations behind the undular bore region 
can be found in 58j, where coUisionless shocks in a two-temperature plasma were studied. 
The trailing edge position of the undular bore in the numerical solution is found by using 
a linear approximation of the amplitude profile, which will be shown in the next section to 
be consistent with the asymptotic modulation solution. 

Although, as we will show in the next section, from an analytical point of view the linear 
dependence of the wave amplitude on position is valid only in the vicinity of the trailing 
edge of the undular bore, the linear numerical interpolation of the undular bore envelope 
turns out to be quite accurate for a significant part of the undular bore (for moderate 
values of the jumps across the bore) and can be used for the determination of the trailing 
edge position in the numerical solution despite the fact that the rear end of the undular 
bore is visually obscured by the above-mentioned linear oscillations (see figure 5). 

The zero-amplitude point in figure 5 corresponds to a "degenerate contact discontinu- 
ity" which moves with the velocity Uc = < . This point (a characteristic in the x — t 
plane) is an analog of the contact discontinuity of classical gas dynamics, although in our 
present dissipationless case there is no discontinuity forming here. 

V. MODULATION SOLUTION IN THE VICINITY OF THE 
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TRAILING EDGE 

Next, we study the asymptotic behaviour of the modulation solution in the vicinity of 
the (given) trailing edge s = s"(A) defined by (jKIH) provided that A < Acr- Since the 
amplitude of the oscillations near the trailing edge of a fully developed undular bore is 
small, we can use the asymptotic modulation system Assuming that the mean 

flow variations in the vicinity of the trailing edge of the undular bore are induced entirely 
by the wave motion we introduce the asymptotic decompositions 

f] = r]o + A%{k) + ... , u = Uo + A\2ik) + ... , (56) 

where rjo, Uq is a a constant flow, and substitute them into the asymptotic modulation 
equations (j^ . ((23) • We assume that (this will be confirmed by the actual solution) 

A^d^A^ < A^d^k < d^A^ . (57) 

Then consistently to leading order in A we get 

+ ^_24+54+84+3 

where Kq = k'^rjl/3. Then, substituting (jHHjl . (jKH|) into the modulation equations (j2EI), 

we obtain the classical modulation equations for a weakly nonlinear wavepacket (see jH], 

Ch.l4) 

^ + c.;(fc)^ + = 0{A'd.A\A'd.k) , (59) 

^ + u',{k)^ + ^2{k)^ = 0{A'd,A\ A'd.k) , (60) 

where u!o{k) = UQ{rio,uo, k) is given by Q and 

~ fc(3/tS + /.6 + 14/t^ + 33^g + 9) 

^^^^^ - " 8^o«:g(l + .g)(.^ + 3.g + 3) < ° 

is an effective weakly nonlinear correction to the frequency ujQ{k), including the effect of 
the induced mean flow. 
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The characteristic velocities of the system (fS^. (jUUj) are readily found, and have the 
form 

X± = uj'o± ^JJ^2A + 0{A^) . (62) 

Thus the modulation system is hyperbolic if co'qCJ2 > 0, which is the classical Lighthill 
criterion for modulational stability of a weakly nonlinear wavepacket (see jH], [21] for 
instance). Since here 

we infer that the weakly nonlinear wave trains of the system are modulationally stable 
for all k provided condition (jKTj) is satisfied. 

Next we construct the similarity solution of the system (jUUI) satisfying the bound- 
ary condition 

x = s-t: A = 0, k = k- , (64) 

where s~(A) is defined by (j^Ujl and k~{A) is found from Eq. ()49p evaluated at f/ = A. 
This similarity solution of equations (jUUj) subject to conditions (jMj) . and an additional 
condition dk/ds < distinguishing the right-propagating undular bore, has the form 

k = k- + ^^is-s-) + Ois-s')\ (65) 

A= — ^ -Js - s') + 0(s - s-)\ (66) 

3^u'^{k-)Co,ik-y 

where u!Q{k) and u;2(/i;) are given by Eqs. (j63|l . (j6T| . The asymptotic solution for fj, u is 
given by ()58|) . We see that the solution ()66|1 in view of ()28|) implies a linear growth in 
amplitude a ~ A ~ (x — x~) near the trailing edge of the undular bore for a given t. One 
can see now that conditions ()57p used in the derivation of the weakly nonlinear reduction 
(jSni), (jHIHl of the modulation system are indeed satisfied for the obtained solution 
so that the whole construction is consistent throughout. 
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VI. FORMATION OF A FINITE- AMPLITUDE REAR WAVE 
FRONT 

The existence of a minimum for the traihng edge speed s~ as a function of the depth ratio 
A across the undular bore (see figure 3) suggests that there is a possibihty of wave-front 
formation at the rear end of the undular bore when A = A^r- Indeed, as s~ from ()39|) has 
the interpretation as a hnear group velocity, there is an analogy with the classical situation 
of the formation of the rapid wave front in linear wave theory (see Ch. 11). In that case 
the existence of a real root k = k* for the equation oo'^ik) = implies that the group ve- 
locity has a minimum (or maximum) at k = k*, so that the wave packet cannot propagate 
with a speed lower (greater) than uj'Q{k*)] consequently a wavefront forms beyond which 
the wave amplitude decays rapidly to zero. From the viewpoint of modulation theory, this 
implies a local linear degeneration of the wave number conservation law 

in the vicinity of a critical point k = k*. Of course, as is well known (see 0) this rapid 
change occurs only in the asymptotic theory, the full solution for the wave front transition 
in linear theory is described by an Airy function, and is characterised by an exponential 
decay in the wave amplitude outside the front. 

A nonlinear analog of this classical wave front formation has been recently noted in 
|UU] for an undular bore solution of the Camassa-Holm (CH) equation, which is a certain 
weakly nonlinear uni-directional approximation of the SG equations (0) (see ^.J. For the 
special family of periodic peakon-type solutions of the CH equation, considered in [HHj . 
the modulation system consists of just two equations, say for k and u, and its expansion 
fan solution exhibits a turning point for the nonlinear characteristics of the modulation 
system, so that the wave cannot propagate back beyond some x = x*. Behind this point 
the solution decays exponentially. The turning point for the nonlinear characteristics is a 
nonlinear analog of the minimum for the linear group velocity. Indeed, the linear dispersion 
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relation for the CH equation allows for real roots for the equation uj^lk) = corresponding 
to the the classical linear wave front formation. 

In the present case of the system the second derivative of the linear dispersion rela- 
tion Q never vanishes (see (jUHj) ). However, in the full modulation theory ujq = ujQ{fi, u, k) 
(see Eq. ((221)) and the condition for wave front formation should be formulated in general 
terms of the linear degeneration of the characteristic field, rather than just the simple van- 
ishing of the second derivatives of the linear dispersion relation. The j-th characteristic 
field Xj of the quasilinear system of conservation laws is called linearly degenerate [33] if 



Vj ■ gradAj = , (68) 

where is the right eigenvector corresponding to the eigenvalue A-,-, i.e. the characteristic 
velocity has an extremum in the direction orthogonal to its "own" characteristic direction. 
For the Whitham system (fT7|) . ()18|) the eigenvalues and eigenvectors are defined in terms 
of an equivalent standard representation pUj) by 

{Xj : det(B - AI) = 0} ; B = A^r^ . (69) 

In the SG modulation system, the condition ()68|1 means that the linear degeneration 
takes place at the point of the minimum for the characteristic velocity A^ in the expansion 
fan similarity solution (jS^)- As the solution implies dXj/ds > 0, the lower bound for the 
characteristic velocity range in the undular bore is determined by the trailing edge speed, 
hence the first appearance of the linear degeneration of the characteristic field A-,- as A grows 
is expected at the trailing edge where the wave amplitude vanishes. Thus the occurrence 
of the linear degeneration point on the undular bore profile is indeed determined by the 
minimum of the curve s^(A), 



dA 

which yields Acr ~ 1-43 (see figure 3). 



ds 

0, (70) 
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It is instructive to also obtain the value of A^r directly from the definition (j68j). We con- 
sider the zero-amplitude reduction of the full modulation system with the built-in simple- 
wave relation u = 2(y^ — 1), which is an exact integral of the ideal shallow- water equations 
fl21|) and which is consistent with the relationship ()37|) between u and r] at the trailing edge 
of the undular bore, 

where V{f]), flo{f], k) are given by Eqs. (jUj), (jS)). Thus for the limiting case under study, 

V(f]) 
B= I _ I , (72) 

dQo{ri,k)/dri dflo^r], k)/dk 

and an elementary calculation yields 

_ dno{f],k) dno{f],k)/df] V 

' dk ' ' V v{f])-dno{f],k)/dkJ ^ ^ 



so that the orthogonality condition ()68p takes the form, 

dkdf] dk^ V{r]) -dno{r],k)/dk ' ^ ' 

Note that (f7^ is a "two-dimensional" analog of the condition lo'^ik^ = for the single 
wavenumber conservation law (jHZj). After substitution of Qo{f],k) and V{f]) Eq. (fTIj) be- 
comes 

- 8a^ - Ua^ + 2a^ + 2a + 2 = (75) 

where a = {1 + k'^ff /'i)^^/'^ < 1. The only real root of ()75|) . which is less than unity, is 
ao ~ 0.661. The intersection of the curve a{fi,k) = with the integral ()49|) relating k 
and f] at the trailing edge of the undular bore yields the critical value of the density ratio 
across the simple undular bore and the corresponding wavenumber at its trailing edge (see 
figure 5): 

Ac^f« 1.430, fc- ^ 1.376 (76) 
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Figure 6: The intersection of the hnear degeneration curve a{ri,k) = ao (sohd hue) with 
the dependence k{f]) at the traihng edge given by (jl^ (dashed hue) yields the critical value 
of A ^ 1.43 

The corresponding values of other parameters at this point are: ^ 0.392, ~ 0.737. 
It is readily shown from the ordinary differential equation for k{fj) that the linear 
degeneration condition (f7^ considered for the trailing edge of the undular bore is equivalent 
to the condition (f7n|) . 

The numerical solutions for different A show that, as A increases further beyond Acr, 
the point of linear degeneration shifts towards the leading edge of the undular bore so that 
a "partial undular bore" with a finite-amplitude rapidly-varying (on the modulation scale) 
rear wave front forms (see figure 7). To determine the speed of a rear finite-amplitude wave 
front in the partial undular bore of this type one should use the general expression (jF)H|l 
rather than its zero-amplitude version ()74|) . An important consequence of the occurrence 
of the rear wave front is that the Gurevich-Pitaevskii type formulation used so far, is not 
applicable for A > Acr as it is based on the assumption of a fully developed undular 
bore with a zero- amplitude trailing edge, which in its turn implies the natural continuity 
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Figure 7: The supercritical, partial undular bore for A = 3.5, t = 120. Initial step data: 
T]- = A, u- = 2{^/A - 1), r/+ = 1, u+ = 

matching conditions for the mean height and velocity at the trailing edge. Instead, 
for A > Acr the conditions at the wave front should presumably be formulated in terms 
of the original Whitham shock conditions for the modulations j7], [S] following from the 
averaged conservation laws (fTTjl . (fT^ . 

Thus, we reach a different problem formulation for A > Acr, and this explains some of 
the discrepancy in figure 3 between the numerical and modulation solutions for the leading 
edge s^(A) (and correspondingly for the amplitude a^(A) in figure 4) as the modulation 
solution constructed in Section V is based on the assumption of a fully developed undular 
bore with the natural boundary conditions satisfied at its edges. Curiously, the numer- 
ical curve for the trailing edge in figure 3 constructed using the formal linear approximation 
for the undular bore envelope, a ~ x — x~, still agrees very well with the formal modulation 
solution for the zero-amplitude trailing edge despite the fact that such a trailing edge itself 
is already nonexistent! Thus, the formal modulation solution (j50|) . (j5T| can still be used 
to estimate the width of moderate-amplitude undular bores. Of course, one should bear in 
mind that this discussion does not essentially affect the shallow-water application of the 
system ((T)), as then undular bores exist only for A ^ 1.3 (PP, jH]), and for such undular 
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bores the modulation theory is in exact correspondence with the full numerical solution. 

Although the detailed description of the structure of supercritical, partial undular bores 
is beyond the scope of the present work, one can infer some qualitative features of such 
undular bores from an asymptotic analysis of the modulations in the vicinity of the trailing 
edge for the undular bore corresponding to A ^ A^r, i.e. in the "nearly-critical" configura- 
tion. For A ^ Acr the rear wavefront is not yet pronounced and the wave amplitude close 
to it is still small as in the vicinity of usual trailing edge. The asymptotic behaviour of 
the modulations near such an emerging wave front, however, will differ from that for the 
trailing edge for "normal" undular bores described in Section V in that the linear degen- 
eration of the characteristic field at the wave front will result in a square root variation 
of amplitude with position, a ~ (x — x^Y^'^, rather than the universal linear behaviour 
a ~ (x — as in (j66|) . This can be understood by considering the characteristic ve- 
locities of the reduced modulation system and ()fi()|l near the turning point where the 
system is close to a linear degeneration, i.e. when co'q(/c) = O^A^) <^ 1. Indeed, in this 
situation the corresponding expansion of the eigenvalues X± (see ()62|) ) would have the form 
A± = i^o(^) + ^(^^)? which implies the square root behaviour of the self-similar solution 
for the amplitude a ~ A ~ (s — One can expect qualitatively the same behaviour 

for the full asymptotic modulation theory described by the system P^ - (P7j) . for which the 
condition for linear degeneration is more complex (see ()74|1 ) . This explains the visually 
noticeable change in the undular bore envelope shape from nearly linear to convex when A 
increases (cf. figures 5 and 7). While for A — A^r = 0{1) one should not expect the same 
square root asymptotic behaviour of the amplitude near the (finite-amplitude) wavefront, 
one can see from figure 7 that the supercritical undular bore retains the same convex en- 
velope shape. We note that a qualitatively similar type of wave envelope behaviour occurs 
for undular bores for the CH equation, in which case the rear wavefront is always present 

As we have already mentioned, the continuous modulation solution does not describe 
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the rapid exponential decay of the wave amphtude at the rear wavefront (see [60 ). One 
can, however, retain the modulation description by introducing appropriate discontinuities 
in the modulation variables using the conservative form of the modulation equations (fTTjl 
and ()18|) . Such formal, discontinuous modulation solutions were proposed by Whitham in 
[7j, but have not been employed so far. 

The presence of jumps for modulations at the rear wave front implies that the Riemann 
invariant condition ()37j) . which is a consequence of the continuous characteristic matching 
|44j . jini is no longer valid. As a result, for A > A^r this condition no longer prevents 
the generation of a backward rarefaction wave. This is clearly seen in figure 7 where 
numerical solution is presented for the decay of a large initial step with A = 3.5 with 
the values r/^ = A and u~ initially satisfying the Riemann invariant condition u~ = 
2{y/T]~ — l). Surprisingly, as our analysis of numerical solutions shows, these "supercritical" 
undular bores apparently follow the classical shock curve (jHHj) . which is quite unexpected 
for dissipationless wave dynamics. Indeed, the classical shock curve occurs under these 
conditions where the transition zone of constant width forms either due to dissipation, 
or due to the existence of a steady smooth bore solution for the original equations (see, 
for instance, the internal bore solution in |SZI)- As our numerical solution shows, the 
supercritical undular bores for the SG system are unsteady, so the mechanism leading to 
the occurrence of the classical bore conditions is not clear. We plan to investigate this 
effect in a separate study. 

In conclusion, one should mention that partial undular bores of the above described 
type are different from the partial undular bores occurring in initial-boundary value prob- 
lems for nonlinear dispersive systems (see for the description of such a partial undular 
bore in KdV theory). Indeed, the latter occur owing to the presence of the boundary condi- 
tions and retain the main qualitative properties of regular undular bores, while the former 
are characterised by the formation of a linearly degenerate wavefront (which is absent in 
the KdV theory) and have a qualitatively different envelope shape. 



39 



VII. DAM-BREAK PROBLEM 

We next consider the classical dam-break problem. Let a separator (a dam) hold back 
water of some depth, say Aq > 1 while the water in front of the dam has unit depth. The 
dam breaks at t = and releases the water downstream. Thus, we are dealing with the 
decay of the initial discontinuity: 

t = : ?7 = Ao > 1, M = for X < 0; r/ = 1, m = for x > (77) 

for the equations (jT}. We first assume that Aq is not too large, so that the generated 
undular bore is "subcritical" , i.e. fully developed. Since the discontinuity (fTTj) does not 
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Figure 8: Numerical solution of the dam-break problem for rj. Initial total depth jump 
Ao = 1.8, t = 150 

satisfy the simple undular bore transition relation (jHTf) . two waves must occur after the 
dam breaking. Simple analysis shows that the relevant combination consists of a right- 
propagating simple undular bore, and a left-propagating centred rarefaction wave (see 
figures 8,9). Then the depth ratio across the undular bore is found from the intersection of 
the transition curve (|37jl and the centred left-propagating rarefaction wave curve v^Aq = 
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Figure 9: Numerical solution of the dam-break problem for u. Initial total depth jump 
Ao = 1.8, t = 150 



u~ /2 + v^r/^, where ?7~, u~ are the parameters of a constant mean flow in the region 
separating the undular bore and the rarefaction wave. As a result we get 

( /7\~ -I- 1 

V- = \ ^ , u- = 2(^-1). (78) 

For instance, if Aq = 1.8 as in figures 8 and 9, then from ()78p we have t]~ = 1.37, = 0.34, 
which agrees with the full numerical solution (of course, this illustration concerns the 
mathematical model itself, and is beyond its scope of applicability to actual shallow-water 
undular bores). 

It is easily seen from Eq. (f75j) that r/^ < Aq. The critical value of Aq corresponding to 
f]' = Acr = 1.43 is 1.94 so that for Aq > 1.94 we get a partial undular bore satisfying the 
jump condition To get the dependence of the lead solitary wave amplitude on the 

initial depth ratio one should set Eq. (fTHj) into the equation For small 5o = Aq — 1 

we have the expansion 

a^ = So- + 0{6l) , (79) 

To compare the small-amplitude expansion ()79|) with the corresponding weakly nonlinear, 
uni-directional KdV result one should use the "effective" initial jump 5 = ?7~ — 1, where ri~ 

41 



is calculated using Eq. (fTHj) . as the value of an initial discontinuity in the KdV equation. 
One can then see that the classical KdV result a+ = 25 agrees with the first term of the 
expansion (ffnjl . 

The distinct oscillations generated at the right boundary of the rarefaction wave in the 
numerical plots in Figs. 8,9 are due to resolution of the weak discontinuity which would be 
present here in the rarefaction wave solution for the dispersionless shallow-water equations. 
These linear oscillations are universally described by the Airy function and have been stud- 
ied in j^HI- Unlike the linear oscillations just behind the undular bore, which decay with 
time as the amplitude of the oscillations generated near the right boundary of the 
rarefaction wave decay as 

VIII CONCLUSIONS 

We have studied the undular bore transition connecting two different constant basic 
states in the fully nonlinear shallow-water model described by the SG equations ((H) derived 
using shallow-water approximation in the full Euler system for irrotational flows. The 
main feature of the undular bore transition is its unsteady character, due to the absence of 
dissipation in the model. Such an expanding undular bore can be represented as a slowly 
modulated periodic wave with the modulations governed by the Whitham modulation 
equations. The wave is confined to an expanding interval s^t < x < s^t, where and 
are the speeds of the leading and the trailing edge respectively. The modulation is such 
that close to the leading edge the wave assumes the form of successive solitary waves while 
close to the trailing edge it degenerates into small-amplitude harmonic wave. 

Using some special properties of the Whitham equations in the small-amplitude (linear) 
and small-wavenumber (solitary wave) limits we have obtained exact analytic expressions 
for the speeds and the lead solitary wave amplitude a"*" as functions of the depth ratio A 
across the undular bore. These results differ from those derived previously for the weakly 
nonlinear case described by the KdV equation. In particular, the fully nonlinear theory 
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predicts a substantially narrower undular bore transition even for moderate values of the 
depth ratios across the bore. It is also shown that one should use the zero-jump condition 
for the ideal shallow-water Riemann invariant across the bore, rather than classical jump 
conditions applicable to an established steady-state undular bores with small dissipation. 
Our analysis is performed for the whole range of A. It is shown that a critical value exists 
Acr ~ 1-43 corresponding to the minimum of the relation s~{A), so that in undular bores 
with A > Acr a rapidly-varying finite-amplitude rear wave front forms instead of a zero- 
amplitude trailing edge, characteristic for small-amplitude undular bores. Such an effect 
is absent in weakly nonlinear KdV theory, but has recently been observed in the more 
sophisticated Camassa-Holm model. 

Parallel to the modulational analysis, we have carried out some direct numerical simula- 
tions of undular bore development in the fully nonlinear shallow- water model (^. Excellent 
agreement of the modulation solution with the parameters drawn from the full numerical 
solution has been demonstrated for "subcritical" (A < Acr) undular bores. The occurrence 
of a rear finite-amplitude wave front for A > A^^ predicted by the modulation analysis is 
also confirmed. Such an agreement can be regarded as a striking confirmation of validity 
of the Whitham modulation theory in unsteady fully nonlinear wave problems, where the 
exact methods of integrable soliton theory are often not applicable. 
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Appendix A: Derivation of the SG system 

For convenience, we here produce a brief summary of the derivation of the system 
originally derived in jHO]. The equations of motion for fully nonlinear irrotational water 
waves are, when expressed in long-wave coordinates X = ex,T = et, z, 

eVxx + 4>zz = 0, 
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w = at z = , 
1 

Vt + (PxVx = -w at z = ri, 
e 

(/)T + -{^/^ + u^^} + ?7 = at z = rj . 
Here the horizontal and vertical velocity are given by 

e 

First we have the exact equation for conservation of mass, which is the first equation in 
©, 

VT + {Ur])x = , 
rv 

where r]U = udz . 







Next we expand in powers of 



F{X,T)-e^-Fxx + 



z^ 

so that u = Fx — e^—Fxxx + 



and w = —ezFxx + ■ • • , 

2 

Also U = Fx - e^^Fxxx + ■ ■ 
o 

Systematic replacement of Fx with the mean flow U yields 

f 
6 

n2 



Fx = U + e'''lrUxx + 



u{z = 7]) = U - e^^Uxx H , 

w{z = rj) = erjUx + ■ ■ ■ ■ 
Hence the Bernoulli condition at the the free surface becomes 

Ft - e^'^UxT + \{U^ - ^^"^UUxx + eWU^} + r/ = . 
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Finally, differentiation with respect to X yields 

Ut + UUx + Vx = e'^ivHUxT + UUxx - U'^)} , 

or] 

This is the second equation in ([T]). Note that only terms of O(e^) are needed, and that w 
is linear in z to this order. 

Appendix B: Numerical Method 

During the development of modulation theory for the fully nonlinear shallow-water 
equations, comparisons have been made between predictions of this theory for various 
properties of an undular bore and the results of full numerical solutions of the equations 
under study. Here, a brief outline will be given of the numerical method used to solve the 
system (P). 

The first of the equations was solved using the Lax-Wendroff method (see [HII for 
instance). To solve the second equation, the time and space derivatives were approximated 
using standard second order, centred finite differences. This resulted in a stable numerical 
scheme that was second order in both the time At and space Ax steps. Due to the {rfuxt)x 
term in the second equation ((T)), solving for u required the solution of a tri-diagonal system 
of equations for u at the new time step. In the numerical solutions shown in the present 
work. Ax varied from 10~^ to 10^^ and At varied from 10^^ to 10~^. 
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Figure Captions 

• Figure 1. 

Oscillatory structure of the undular bore evolving from an initial step (dashed line). 

• Figure2. 

Riemann invariant (dashed line) and the shock-wave (dotted line) transition curves. 

• Figure 3. 

The leading s"*" (upper curve) and the trailing s~ (lower curve) edge speeds vs the 
depth ratio A across the simple undular bore. Dashed line: modulation solution ()50p. 
(15111 : Diamonds: values of s^, s^ extracted from the full numerical solution; Dotted 
line: the KdV modulation solution 
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• Figure 4. 

The lead solitary wave amplitude vs depth ratio A. Solid line: numerical solution; 
dashed line: modulation solution (jSH). 

• Figure 5. 

The numerical solution of the system for the initial depth ratio A = 1.2, t = 250. 
Initial step data: r]^ = A = 1.2, u~ = 2{y/r]^ — 1), 77+ = 1, u+ = 0. 
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